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Abstract 

A theory for the generation of current in a toroidal plasma by radio-frequency waves is presented. 
The effect of an opposing electric field is included, allowing the case of time varying currents to be 
studied. The key quantities that characterize this regime are identified and numerically calculated. Circuit 
equations suitable for use in ray-tracing and transport codes are given. 
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I. INTRODUCTION 



In recent years there has been considerable interest in generating steady-state currents in a plasma 
with rf waves. In particular, it was predicted^ that these currents could be efficiently generated by waves 
whose phase velocities are several times the electron thermal speed. This prediction has been confirmed 
by numerous experiments in which the current was driven by lower-hybrid waves. These results allow 
us to contemplate a steady-state tokamak reactor in which the toroidal current is driven by lower-hybrid 
waves. This is an attractive proposition not only because of the advantages inherent in steady-state 
operation (less thermal stress, higher duty cycle, etc.), but also because it opens up the possibility that 
the ohmic winding of the tokamak can be eliminated entirely, leading to a cheaper and more compact 
reactor. This latter possibility can be realized if rf waves are successful not only in sustaining the plasma 
current "steady-state current drive," but also in increasing the plasma current "rf current ramp-up." In fact, 
experiments have demonstrated that this too is possible. From a theoretical point of view, the important 
additional ingredient in these experiments is the dc electric field which opposes the increase of the plasma 
current. The electric field is also present in schemes where the rf is used to recharge the transformer at 
constant current. 

Recently, we presented a theory for rf current drive in the presence of an electric field.^ This theory 
predicted that rf energy could be efficiently converted to poloidal field energy if the wave phase velocity 
were approximately equal to the electron runaway velocity. This theory has been compared"^ with data 
from the PLT experiment,^ and excellent agreement is found. In Ref . 2, the linearized Boltzmann equation 
was approximately solved by integrating the corresponding Langevin equations using a Monte-Carlo 
method. In this paper, we use a more elegant theory to calculate the efficiency of the current ramp-up 
based on an adjoint formulation for the Boltzmann equation.^ Although the Umits of validity of this 
theory are the same as for Ref. 2, this theory is more amenable to accurate evaluation on a computer, and 
it is more easily extended to include effects which are omitted here. 

Let us begin by reiterating the physical picture given in Ref. 2. Consider an electron travelling in 
the positive direction at several times the thermal speed and which has just absorbed an incremental 
amount of rf energy. Suppose there is an electric field tending to decelerate this electron. The question is: 
Where does this incremental energy end up? If the electron is slow compared to the runaway velocity, the 
electron slows down primarily due to collisions and so the rf energy goes to bulk heating. On the other 
hand, if the electron is fast, the electron is slowed down by the electric field. In this case, the rf energy 
is coupled to the plasma circuit and appears as poloidal field energy. Unfortunately, fast electrons have a 
high probability of pitch-angle scattering into the reverse direction and running away. A runaway electron 
drains energy out of the electric field, leading to a degradation of the ramp-up efficiency. However, there 
is a window around the runaway velocity where the electrons are slowed down principally by the electric 
field and yet where the probability of running away is very small. This is the favorable regime in which 
rf energy can be efficiently converted to poloidal field energy. 

From the foregoing discussion, we see that two ingredients are needed for an accurate theoretical 
treatment of this problem. First, the electric field must be treated as large. In the efficient regime, the 
force on the electron due to the electric field must be comparable to that due to collisions. Second, a 
two-dimensional treatment is required. Analyses based on a one-dimensional Fokker-Planck equation 
do not predict the important physical phenomenon of rf-generated reverse runaways. 

We briefly review the history of theoretical studies of current drive in the presence of an electric field. 
The earliest theoretical studies^ of lower-hybrid current drive assumed that there was no electric field 
in the plasma. This is the appropriate limit for a steady-state reactor. However, some of the early ex- 
periments conducted to verify the predictions of the theory were conducted in regimes where the ohmic 
electric field was still present. This prompted a series of papers^"^" deahng with rf current drive in the 
presence of an assisting electric field (i.e., the electric field and the rf both accelerate the electrons in 
the same direction). The principal focus of these papers was the calculation of an enhanced runaway 
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rate when the phase velocity of the waves is in the neighborhood of the runaway velocity. An opposing 
electric field was treated by Borrass and Nocentini^^ within the framework of a one-dimensional Fokker- 
Planck analysis. A similar model has been employed more recently by Liu et al}^ As noted above, such 
an analysis cannot include rf-generated reverse runaways. Also due to the crudity of the one-dimensional 
equation, the results are only accurate to within a factor of order unity even when no runaways are gen- 
erated. A two-dimensional treatment of the problem in the small electric field limit has been given by 
Start^^ for the case of current drive by electron-cyclotron waves. This work neglected electron-electron 
collisions. This defect was removed and the results generaUzed to arbitrary current-drive methods by 
Fisch^'' employing an adjoint formulation. This work yields accurate results when the electric field is 
small. However, the results are inapplicable in the regime of efficient ramp-up where the phase velocity 
is comparable to the runaway velocity. Our earlier paper^ was the first to combine a two-dimensional 
treatment with a large electric field. This paper allowed an accurate calculation of the ramp-up efficiency 
in cases of practical interest, and identified the regime in which high efficiencies can be expected. The 
present work is a continuation and expansion of that earlier paper. Besides using the more sophisti- 
cated method'' for solving the Boltzmann equation, we endeavor to give the results in a form that allows 
both easy comparison with experiments and easy implementation within the framework of ray-tracing or 
transport codes. 

The paper is organized as follows: We begin with the linearized Boltzmann equation for the perturbed 
electron distribution in the presence of an electric field and an rf source (Sec. II). Some approximations 
and normalizations are made to reduce this equation to a more manageable form. The use of the adjoint 
method^ for solving the resulting equation is described (Sec. III). Next (Sec. IV) the adjoint equation is 
solved numerically to give the runaway probability and the Green's function for the current. The latter 
quantity is reduced to a simple form which involves just two functions of velocity. An expression for 
the total current density is given (Sec. V) and this is put into a form which is easy to calculate. How the 
rf-driven current interacts with the electric field to produce poloidal field energy is considered (Sec. VI), 
and the results are applied to experiments. 



II. BOLTZMANN EQUATION 

Consider a uniform electron-ion plasma, initially at equiUbrium. For f > 0, it is subject to an electric 
field E{t) and a wave-induced flux S(v, t). We will take the ions to be infinitely massive, so that they 
form a stationary background off which the electrons colhde. If the electric field and the wave-induced 
flux are weak enough, the electron distribution remains close to a MaxwelUan for £ <T where £ is the 
energy of an electron ^mAp . Substituting f = fm + fi into the Boltzmann equation for the electron 
distribution / and linearizing then gives 

where 

( 772 \ ^/^ 
■^^j expi-S/T), 

and 

c{f) = c(/, /„o + c(/„, /) + c(/, fi) 

is the linearized collision operator. Here q, rn, n, and T are the electron charge, mass, number density, 
and temperature. Note that q carries the sign of the electron charge (i.e., q = — e). 

This equation is to be solved with initial condition /i(v, t — 0) — 0. We demand that the subsequent 
evolution of /i be such that it be orthogonal to 1 and £, i.e., that it have zero density and energy. The 
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zero-density condition is satisfied with fi = since all the terms in Eq. (1) are particle conserving. The 
zero-energy condition gives an equation for the time evolution of T 



3 dT f [' 
n— = / mS • V d^v + E • / gv/i (fv. 
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The two terms on the right-hand side represent the heating due to the waves and due to E • J. 

We now make three simplifying assumptions: we assume that fi is azimuthally symmetric about 
the ambient magnetic field; we take the electric field to be constant and in the direction parallel to the 
magnetic field E = _Bv|| ; and we restrict our attention to those cases where S is only finite where v ^ Vt, 
where = T/m is the thermal velocity. We may then solve Eq. (1) using the high-velocity form for C: 



c(/) = r 



1 d ^ ^ l + Z d ^^)^ f 
dv 2v^ djjL 



where ^ = v\\/v,T = nq* In A/47reQm^, eg is the dielectric constant of free space, In A is the Coulomb 
logarithm, and Z is the effective ion charge state. We have included pitch-angle scattering and frictional 
slowing down, but ignored energy diffusion. In the problem of steady-state current drive, the energy 
diffusion term introduces corrections of order (vt/v)'^. Another term neglected in this approximate col- 
lision operator is the effect of the Maxwellian colliding off the perturbed distribution, C(f„i, /). The 
corrections due to this term^^ are of order {vt/vY . With these approximations, the collision operator 
does not depend on the electron temperature T. Formally, we may derive the form for C by taking 
T ^ 0. In this limit, we have n5(y). 

It is convenient to introduce some normalizations. The runaway velocity Vr is that velocity at which 
collisional frictional force equals the acceleration due to the electric field 

Vr = -sign{qE)^JmT/\qE\. 

Notice that the sign of Vr is opposite to the direction in which electrons run away. The Dreicer velocity^^ 
is given by —\/2 + Zvr- Similarly, we define a runaway collision frequency 

Vj. = r/ \Vr\^ . 

The normalized time and velocity are given by 

T = Vrt 

and 

U = v/Vr- 

The components of v have to be normalized with care: u± = v±/ \vr\, M|| = v^\/vr, and u = v/ \vr\. 
This implies that = sign(u,.)?;||/?,' so that the conversion of the pitch-angle variable /i in v space 
to that in u space involves multiplication by sign(t!r). Other quantities are normahzed in a similar 
way; however, we shall use the same symbols as for the urmormalized quantities. Thus, the distribution 
functions /i and /,„ are normalized to n/ \vr\^, the rf-induced flux S to nVrVr/ \vrf, etc. Under this 
normalization, Eq. (1) then becomes 

lf,+Dif.) = -l.S, (2) 

with /i (u, T = 0) = and with the operator D defined by 

d Id 1 + Z d ^ d_ 
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Equation (2) is singular at the origin. However, because /i in Eq. (2) describes a physical particle 
distribution, it obeys a particle conservation law near the origin. We therefore require that, close to 
u = 0, /i(u) = N{t)6{u) with N{0) = and 



— — = lim / 2nti[u,fj.)aiJ.. 



Equation (2) depends only on a single parameter Z. The dependence on the electric field E can be 
normaUzed away, since the electric field defines the only natural velocity scale in the problem, Vr. 

Equation (2) is amenable to various methods of solution, and it is instructive to review these before 
describing the method used here. The most straightforward approach is to integrate Eq. (2) directly on 
a computer. This method allows S to be determined directly in terms of the electron distribution /. 
However, a thorough understanding of the problem requires that many different forms of S be used. 
Therefore, this procedure is costly because the several parameters used to specify S must be scanned. 
This is essentially the method adopted in the early numerical studies of steady-state current drive by 
lower-hybrid waves. 

The situation is improved to some extent by noting that Eq. (2) is a linear equation for /i. It may be 
solved in terms of a Green's function g given by the equation 

(J;^ + D^g{u,T;u')=0 (3) 

with g{u, T = 0; u') = 5{u — u'). The electron distribution is then given by the convolution 

h{n,T)= j^T j (iVS(u',r')-^5(u,T-T';u'). 

This approach reduces the problem to the determination of a single function g of two vector arguments 

(u and u') and a one scalar argument (t). However this is still a daunting computational task. 
A closely related technique is to formulate the problem as a set of Langevin equations, 

du 1 
dr u-' 



dr u 



(4) 



where the pitch-angle scattering is represented by the stochastic term A{t). Assuming that u{t) =u and 
Ii{t) = /U are given (i.e., non- stochastic), then A{t) satisfies 

(5) 

(A(r)A(r')) = ^(l-/xV(T-TO, 

where the angle brackets denote averaging over the ensemble defined by all the realizations of A. Con- 
sider following a particular electron using Eqs. (4). Suppose that the electron is observed to travel with 
velocity u' at t = 0. Then ^(u, r; u') d^n is the probabiUty that the velocity of the electron at time r is 
in the volume element rf^u located at u. In Appendix A, it is shown that this conditional probability g 
satisfies Eq. (3). Thus the solution to Eq. (3) can be found by determining the distribution of a large num- 
ber of electrons obeying Eqs. (4) with initial conditions with u(t = 0) = u'. Consequently, moments of 
g can be determined by ensemble averages of the Langevin variables. For example, the current given by 
Eq. (3) may be found by 

d^uung{u,n,T;u' = {u{t)h{t)) . 



I- 
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Equations (4) may be integrated numerically by noting that 

/■t+At 



/: 



A{t') dr' 



should be picked from an ensemble with mean — (1 + Z)/xAr/u^ and variance {1 + Z){1 — /x^)At/m^ 
where u and /i are the values of those variables at time r. As long as Ar is sufficiently small, further 
details about the distribution of A are unimportant. 

Now Eqs. (4) are the equations solved in our earher paper.^ This shows the exact equivalence between 
the approach adopted there and that employed in the present work. Because the Langevin equations 
describe the electron behavior in a shghtly more physical manner, they often help in the interpretation 
of the solutions to the Boltzmann equation. This is especially true when some electrons run away. The 
Langevin equations are also very easy to solve numerically by a Monte-Carlo method (as was done in 
Ref. 2), although their solution tends to be much more costly than just solving Eq. (3) directly. 

Equations (4) however may be easily solved analytically in the limit u ^ (this is equivalent to 
taking the Umit — > 0). Taking an ensemble average of the equations, we obtain 

du 1 



dT 



U 



2 ' 



__i + z 

dT ~ «3 W- 

Note that u is not a stochastic variable in this limit. Consequently, the hierarchy of moment equations 
may be closed at this point. These are the slowing-down equations solved by Fisch and Boozer^^ to give 
the current moment of electron distribution (u/x) = u{ii). This shows that the approach used in that 
paper is equivalent to solving the Boltzmann equation. 



III. ADJOINT METHOD 

The methods for solving the Boltzmann equation described in the previous section all entail a large 
amount of computation. (An exception is the limit E ^ Q, when the ensemble-averaged Langevin 
equations can be solved analytically.^^) The problem with these methods is that they are all capable 
of giving the electron distribution function /i. Since, in many cases, we are only interested in specific 
moments of /i, we may hope to reduce the computational requirements substantially by using a method 
that gives only those specific moments. Suppose we wish to determine a particular moment of /i , namely 



H{t) = J d3u/lo(u)/i(u,T). 



(For instance, the current density would be given by ho{u) = Let us define the corresponding 
moment of the Green's function 

h{u',T)= J d''uho{u)g{u,T;u'). (6) 
The moment H is then given in terms of h by 

H{t) = ^ dT' J d^u S(u, r') • ^h{u, t-t'). (7) 

What is needed is some method of calculating h which doesn't involve having to find g. This is 
provided by the adjoint formulation of Fisch. ^ He shows that h{\i, r) satisfies 

^h + D*{h) = 0, (8) 
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with h{u, r = 0) = ho{u) and with the operator D* defined by 

D* = — + -—- — (1 - m')— 
du\\ vP' du 2u^ dfx dji 

The singularity at the origin is handled by the boundary condition /i(u = 0, r) = 0. The operators D and 
D* are adjoint operators, so that 

j {hD{f)-fD*{h)) (i3u = 

for all /(u) and h{u) satisifying /(u — > oo) = and h{u = 0) = 0. 

Similar techniques were introduced earlier by Antonsen and Chu^" and by Taguchi^'^ for the study of 
steady-state current drive. The significant improvements afforded by Ref. 5 are the abihty to determine 
arbitrary moments of /i and the inclusion of the time-dependence of /i. Both of these are important in 
the problem of current ramp-up. 

From the relation between the two Green's functions h and g, we see that h has a simple interpretation. 
Equations (3) and (4) describe the evolution of a group of electrons released at r = at velocity u'. Let 
us suppose that we are interested in the current density so that /io(u) = Then h{u' , r) gives the 
mean current carried by those electrons a time t later. How Eq. (8) works is easily seen by taking 1 
so that the electron only experiences the electric field. In the Boltzmann equation, the electrons have 
slowed down to u' — ru|| at time r. Correspondingly in the adjoint equation, the initial condition ho is 
transported in the reverse direction so that h{u' , r) = ho{u' — tU||). Thus at time r, we are provided 
with information about the electrons in their current location. 

Solving the Boltzmann equation by means of the adjoint formulation results in a great simplification 
of the problem. The adjoint equation (8) is an equation of equal complexity to the original Boltzmann 
equation (2). However, by solving Eq. (8) for a particular initial condition, we can find the corresponding 
moment of /i using Eq. (7) for any driving term S. 

The proof that h is given by Eq. (8) is most easily carried out by assuming that Eq. (8) holds and 
then by proving that h is related to the general Green's function g by Eq. (6). Consider the equation for 
g{u,T';u'), 

^+D^g{u,r';u') = 0. 

We multiply this equation by h{u, r — r'), integrate over all velocity space, and use the adjoint relation 
between D and D* to give 

J d''uh{xi,T - T')-^g{u,T';u') + g{u,T';u')D*{h{u,T - t')) = 0. 

Substituting from Eq. (8) and integrating in r' from to r, gives 
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d^uh{u,T- T')g{u,T';u') 



= 0. 

r'=0 



If we evaluate this expression using the initial conditions for g and h, we obtain Eq. (6). 

Up until now, we have assumed that all the equations are solved in an infinite velocity domain. 
This is not a convenient formulation for numerical implementation where, necessarily, we wish to solve 
equations in a finite domain. Here we shall only solve Eq. (8) in a spherical domain V such that u < Ub. 
We will choose Ub to be sufficiently large that the interesting physics where the electric field competes 
with the collisions happens inside V. Outside V collisions may be ignored and the electrons are merely 
freely accelerated by the electric field. We must impose boundary conditions on S, the boundary of V. 
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Again we follow the treatment given in Ref. 5. We begin by noting that both Eq. (2) and Eq. (8) are 
hyperbolic in the u direction. The boundary therefore divides into two pieces depending on whether 
the characteristics enter or leave the domain. We define Sin (resp. Sout) as that portion of S on which 
— /X — Xjv^ < (resp. > 0). The characteristics of Eq. (2) enter V on Sjn and leave on Sout. while 
those of Eq. (8) enter V on Sout and leave on Sin. Boundary conditions must be specified where the 
characteristics enter the domain V. For «(, sufficiently large, the solution beyond S may be determined 
by ignoring collisions. Thus 

/i(u,r)=/i(u + ru||,0)=0 

for u on Sin and 

/i(u,r) = h{u — TU||, 0) = ho{u — ruy) 

for u on Sout- (If energy scattering had been included in the collision operator, the equations would 
revert to parabolic, and boundary conditions would have to be specified over the whole of S. However, 
there would be a boundary layer where the characteristics of the approximate equations are outgoing, and 
the boundary conditions here would only weakly affect the solution in the interior of V.) 

Although Eq. (8) was derived under the simplifying assumptions that the electric field was constant 
and the high- velocity form of collision operator is vaUd, the adjoint method as described in Ref. 5 appUes 
equally well without such restrictive assumptions. Thus, the equation adjoint to Eq. (1) reads^ 

9 _sEL:n.^_ c*)hi.,t'-,t) = + s, (9) 

of m ov J 

where C* is the operator adjoint to C. Since the full linearized collision operator is self-adjoint, we have 
C*(h) = C{fmh)/fm- Equation (9) is to be solved with the initial condition /i(v, t' = 0;t) = ho{v). 
We restrict f^h to being orthogonal to 1 and £ and qi and q2 are chosen to ensure that this condition on 
fmh remains satisified given that it is satisfied initially. It can then be shown^ that 



/ 



d3v/i(v,i)/io(v) = j^dt' j <fvS,{w,t') ■ -^h{M,t-t';t), (10) 



where S,(v, t) = S(v, t) + [qEj{t)/rn)j f^. This equation will enable us to incorporate the effects of a 
slowly varying electric field into our analysis. It also makes explicit the additive nature of those effects 
due to the electric field alone (i.e., with S* = (gE/m)/™) and those effects due to the rf combined with 
the electric field (i.e., with = S). Of course, the effects due to the electric field alone are well studied 
and give rise to phenomena such as the Spitzer-Harm conductivity^^ and runaways. 



rV. SOLUTIONS TO THE ADJOINT EQUATION 

Moments of the electron distribution /i can now be calculated by solving Eq. (8) with the correspond- 
ing initial and boundary conditions. In practice, this procedure still offers us too much information. Both 
for a deeper understanding of the underlying physics and for easy implementation in numerical codes, 
the trick is to discover the few important functions by which the major effects can be described. In this 
section, we determine those functions needed for an accurate treatment of rf current ramp-up. 

Let us suppose that rf flux is present only for some finite time. Electrons obeying Eq. (2) then 
eventually suffer one of two fates. Either they run away under the influence of the electric field u\\ 
— oo, or else they collapse into the electron bulk u ^ i). We classify these two groups of electrons as 
"runaway" (subscript r) and "stopped" (subscript s) respectively. In a real plasma, i.e., T ^ 0, even the 
bulk particles will eventually run away. However the time tr it takes for these bulk electrons to run away 
is exponentially large, i.e., log tr ~ {vr/vt)"^. Our analysis is valid for times short compared with the 
bulk runaway time. 
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Runaways are very important in the calculation of the ramp-up efficiency because runaways gain 
energy at the expense of the poloidal magnetic field. Unless they are lost, even a small number of 
runaways can greatly reduce the ramp-up efficiency. Runaways may be defined as those particles with 
u > Wo for T ^ oo where uq is some arbitrary positive speed. (The number of runaways is independent 
of uq.) Therefore their number is given by Eq. (7) with r — > oo and ho{u) = 1 for u > uq and 
otherwise. The Green's function for the runaway number is given by i?(u) = h{u, r — > oo) where R 
obeys 

D*{R{u))=0, (11) 

with boundary condition -R(u) = 1 on Sout- This function is the "runaway probability," the probability 
that an electron initially at u runs away under the combined influence of the electric field and collisions. 

Equation (11) was solved numerically with the boundary at = 10. A term dR/ dr was included on 
the left-hand side, and the resulting equation was integrated until r = 100. A spherical (w, 9 = arccos jj.) 
grid was used with a mesh size of 500 x 100. The equation was integrated with an ADI (alternating 
direction implicit) scheme with a time step At = 0.01. The same method was used to solve the other 
equations given below. 

In Fig. 1, we plot i?(u) for Z = 1. For u < 1, R is identically zero because the magnitude of the 
electrical force is less than that of the frictional force. One of the most important applications of these 
results is to drive current by lower-hybrid waves. In this case S is in the parallel direction and is localized 
near u± = 0. Therefore, we need only know || , u± = 0) which is plotted in Fig. 2 for Z = 1,2, 5, 
and 10. From this plot, we see that we can effectively avoid the creation of runaways by operating with 
waves whose phase velocities he in the range < un < 1.5. 

The next important quantity to determine is the current density carried by /i . This is given (in units 
of qnvr) by Eq. (7) with /io(u) = The Green's function for the current j'(u, r) is therefore given 
by dj/dr + D*{j) = with initial condition j(r = 0) — U|| and boundary condition j — tt|| — t on 
Sout- This is the mean current (in units of qvr) carried by an electron initially at velocity u. In Fig. 3(a), 
we plot j(u, r) as a function of r for u = 5u|| and Z = 1. Because the presence of runaways leads to 
a secular behavior (j ~ r) for large times, it is helpful to distinguish the current carried by stopped and 
runaway electrons. We write 

i(u,r) = (1 - R{u))Mu,t) + J?(u)>(u,r). 

The quantity js (resp. jr) is the mean current carried by an electron given that it eventually stops 
(resp. runs away). An electron at velocity u runs away with probability i?(u). Thus it contributes 
(l — i?(u))u|| to the stopped current and i?(u)M|[ to the runaway current. These quantities are therefore 
the initial conditions to the adjoint equations for (1 — R)js and Rjr respectively, so that 

with js{T = 0) = jr{T = 0) = W|| and js = jr = - T on Sout- 

The stopped and runaway currents js and jV are plotted in Figs. 3(b) and (c) for same case as Fig. 3(a). 
Evidently, j, vanishes for r ^ oo (since the electrons cease to carry any current once they are stopped). 
The time it takes for the electrons to be stopped is of the order of u. Assuming that this time is short 
compared to the time scale for the variation of the rf flux S, we may replace js(u, r) by Ws{u)5{t) 
where Ws(u) = Jg°° js(u, r) dr. The equation for Wg is obtained by integrating Eq. (8) over time to 
give 

D*[{l-R{u))Ws{u)] = (l-i?(u))«||, (12) 

with (1 — R)Ws = on Sout- Ws can be interpreted as the energy (in units of mw^) imparted to the 
electric field by an electron as it slows down. In Fig. 4(a), we plot Ws(u) for Z = 1. In the limit 



9 



U|| oo, collisions are extremely weak, and all of the kinetic energy of the stopped particles goes into 
the electric field, i.e., 

Ws{u\\ ^ oo,u± = 0) — i- 

In the limit u <^ 1, the electric field weakly perturbs the electron motion. Then, Ws is given by the 
theory of steady-state current drive,^^'^° and corrections linear in the electric field are given by the hot 
conductivity.^^ In our notation, these results may be summarized by 

(2 + Z + 3Ai2)u6 

This function is plotted in Fig. 4(b). (An approximation to Wg correct to order is given in Ap- 
pendix B.) This linearized theory, however, is inapplicable for m ~ 1, and the behavior of this function is 
completely wrong for m > 1. 

Let us now turn to the contributions of the runaways to the current. The leading order contribution to 
jr is — T. Let us therefore write 

>(u, r) = -r + >o(u) + j;(u, r), (14) 

where j^(r — » cxd) ~ 1/r and jro may be interpreted as the effective starting velocity for the runaways; 
see Fig. 3(c). The function jVo is given by 

D*{R{u)jro{u))=R{u), 

with boundary condition jVo = W|| on Eout- This function is shown in Fig. 5. For w » 1, the runaway 
electrons are only weakly perturbed by collisions so that jro(u) ~ M ||. Close to u = 1, collisions hold 
back the runaway electrons and jro(u) becomes large. However, it is not very important to know jV-o and 
j'^ very accurately since they are usually dominated by the first term in Eq. (14). We will approximate 
jro(u) by u\\ and will ignore i^(u, t) to give .?V(u, t) ~ M|| — r. 

Finally, we can write the Green's function for the current in an expedient form as 

j(u,r) = (1 - it!(u))W,(u)5(T) +i?(u)(«|| - r). (15) 

In this form, it depends only on two scalar functions of u, namely R and Ws. Approximate fits to these 
functions are given in Appendix B. An easy but important generalization is possible here and that is 
to allow a loss mechanism for runaways. This is done by modifying suitably the term M|| — r. For 
example, if the loss of runaways can be characterized by a loss time rioss. this term should be multiplied 

by exp(-T/rioss)- 

Using this formulation, many other moments of /i may be found. For example, we may wish to 

know the mean perpendicular energy of the runaway particles £±r (in units of mw^) as they leave the 
integration region V. (The loss rate for runaways may depend on this quantity.) This is given by 

D*{R{u)£^r{u)) =0, 

with £±r = on Sout- (This result depends logarithmically on the value of Ub.) We have plotted 
this in Fig. 6. For electrons with «|| > 1 and u± = 0, £±r is about 3. This reflects the necessity for the 
electrons to suffer appreciable pitch-angle scattering if they are to run away. 

V. CIRCUIT EQUATIONS 

When rf energy is injected into a tokamak, it induces a flux S of electrons in velocity space. The 
power deposited per unit volume is then given by 

p,f(t) = j d\mS{w,t) ■ V. (16) 
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Here, prf , S, and the other intensive physical quantities introduced in this section also depend on position 
r. For brevity, this dependence is not shown in the arguments to these quantities. This equation may 
be used in two ways. In detailed studies of rf current ramp-up, based for instance on a ray-tracing 
model, we can estimate S{v,t) on each flux surface by solving a one- or two-dimensional Fokker- 
Planck equation. Equation (16) then gives us the power deposition, p^f. Alternatively, we can take the 
experimental measurements together with an energy balance of the rf energy to give us an estimate of 
Prf . This, together with an approximate knowledge of where in velocity space the rf flux is localized, 
allows us to determine S. In addition to causing power absorption, the flux S leads to numerous other 
effects, such as rf-driven currents, rf-enhanced particle transport, etc. Here our primary concern is with 
the rf-driven current. From Eq. (10), we see that this enters additively to the ohmic current so that the 
total current density is given by the constitutive relation 

J{t) = a{t)E{t) + Jrf(t), (17) 

where a{t) is the Spitzer-Harm conductivity^^ for a Maxwellian plasma characterized by the background 
electron temperature T{t), and Jj-f is the rf-driven current density. Here we have assumed that \vr \ 3> vt 
so that in the absence of any rf we can ignore runaways. Incorporation of this effect merely requires the 
addition of the current carried by the Dreicer runaway electrons in Eq. (17). 

The rf-driven current density is given by Eq. (10) with the h replaced by the current Green's function 
j and with S* = S. Let us begin by writing j in unnormalized units. The form for j given in Eq. (15) 
will be sufficiently accurate for our purposes. Multiplying by qvr gives 

j(v,t) = ^(1 - R{u))Ws{u)S{t)+qR{u)(v\\ + 

where u = v/vr. Here j is now a dimensional quantity, but Wg and R remain dimensionless functions of 
a dimensionless argument. In deriving this form for j we assumed that E and n were constant. We now 
relax this constraint, allowing them both to vary on a time scale long compared to the runaway collision 
time (RecaU that S is also aUowed to vary on the same time scale.) We can then write j as 

j(v, t - t'; t) = (1 - i?(u')) V^.(u')5(t - t') + qR{n') (v\^ +^ f E{s) ds\ , (18) 

where u' = v/t>r(t')- The additional parametric argument t here has the same meaning as in Eq. (9). In 
this form j(v, t — t'\t)\s, the mean current carried by an electron at time t given that it was traveling at 
velocity v at time t' . From Eq. (10), the rf generated current density may now be written as 

JrfW = J d''vS{v,t') ■ ^j{v,t-t';t). (19) 

In order to write in a more useful form, we first define a runaway density rir (in electrons per unit 
volume). This is given by 

9nr{t)_ 1 /•rf3^s(v,i).#i?(u), (20) 



Vr{t) J 



dt Vr{t) J du 
with initial condition nr{t = 0) = 0. Substituting Eq. (18) into Eq. (19), we obtain 

Mt) = Js{t) + Jr{t), (21a) 

where 

Mt) = ^/ rf'vS(v,i) • ^(1 - R{n))W,{u), (21b) 

= ^E{t)nr{t)+q J rf3vS(v,i) • ^R{u)u^^, (21c) 
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with Jr{t = 0) = 0. In Eqs. (20) and (21), u is normalized in terms of the runaway velocity at time 
t, u = v/vr{t). These equations allow the current to be calculated by characterizing the runaway 
population with just two state variables rir and Jr. Equations (16), (17), (20), and (21) suffice to give a 
detailed description of rf current ramp-up. In this form, Eq. (17) is suitable for substituting into a transport 
or ray-tracing code. Furthermore, it would be easy to modify Eq. (20) to include a loss mechanism for the 
runaways. Relativistic effects on the runaways could be included in an approximate fashion by limiting 
\Jr{t)/qnr{t)\ to c the speed of light. Such effects could be treated in a more systematic manner by 
modifying the term in large parentheses in Eq. (18) to read W|| (t) where 




7=^l+p|(t)/m2c2. 

The resulting expression for j(v, t — t';t) is valid for <C and v-^ <C c^. Unfortunately, this is a 
significantly more cumbersome expression from which to calculate Jrf because, in order to determine the 
state of the plasma at a particular instant, the entire runaway distribution must be given (instead of just 

rir and Jr). 



VI. APPLICATIONS 

The circuit equations written in Sec. V allow us to explore how Jrf interacts with the electric field 
to yield an efficient conversion of rf energy into poloidal magnetic field energy. It is helpful to convert 
to extensive physical quantities by assuming that the plasma current is carried in a channel of area A in 
which the plasma properties are approximately uniform. Thus, the total current is given by / = AJ, the 
total rf power deposited in the electrons by Pin = 2'jTRoApj.f (where Rq is the tokamak major radius), 
the loop voltage by ^ = 2'jtRoE, etc. The plasma current is again written as the sum of ohmic and rf 
contributions 

I = V/Rsp + Ir{, (22) 

where Rsp = 2-kRq/A(j is the plasma (Spitzer-Harm) resistance. Faraday's law relates the rate of change 
of the current to the voltage 

V = -Li + V,^t, (23) 

where L is the total plasma inductance, which for simphcity we shall take to be constant, 14xt is the 
voltage induced by the external coils (usually a combination of the ohmic windings and the vertical field 
coils), and I = dl/ dt. Multiplying this equation by I and substituting for I from Eq. (22) gives 

= Pext + Pel (24) 
-n-Sp 

where W = \LP is the poloidal field energy, Pext = Vc^tl is the power coupled from the external 
circuits, and Pei = — V/rf is the power coupled from the rf source into electromagnetic energy. This 
equation describes the energy balance for the poloidal magnetic field. The practical measure of the 
efficiency of current ramp-up is 

W - Pext _ Pel - l^yPsp 

Prf ~ Prf ' ^ ^ 

where Prf is the total rf power injected into the plasma. The rf power absorbed by the electrons Pin is 
related to Prf by Pin = JjPrf where r} is the absorption factor. The determination of r] is beyond the scope 
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of this paper; presumably it can be found by ray-tracing theories or by a power balance. The overall 
picture of the flow of power in an experiment is as follows: Rf power Prf is injected into the machine. 
Of this a fraction rj is absorbed by the resonant electrons; the rest may be absorbed by the ions or by the 
vacuum vessel. A fraction Pei/Pin of this power is then converted into electromagnetic energy. Pext acts 
as another source of poloidal field energy, while the ohmic dissipation V'^/Rsp acts as a drain. From 
this discussion, we see that Pd/Pin describes the "ideal" efficiency of rf current ramp-up. The practical 
efficiency is expressible in terms of this efficiency, rj, and V^/Rsp. 

The determination of Pel/ Pin from Eq. (21) is complicated by the presence of runaways. Runaways 
are deleterious to the ramp-up efficiency since their current is in the same direction as E and so they 
subtract from Pel- For efficient current ramp-up we must either avoid creating runaways by making sure 
S is localized in that region of velocity space where the runaway probabiUty R is small (see Figs. 1 and 
2), or else take steps to lose the runaways. We can approximately treat these cases by taking P = in 
Eq. (21) to give 

Pel _ Jd^uS-dWJdu 
~ /rf3uS-u ■ 

Since this involves the ratio of two integrals over u, the result is insensitive to the detailed form of S. In 
cases of practical interest, we may assume that S is localized in u. Then, we have 

Pel s • dwjdu 
w ^ — 1 ' 

^in S • U 

where u is the normalized velocity of the resonant electrons. 

For lower-hybrid waves we have S = uy and the waves interact with particles through the Landau 
resonance w — fcy = 0, where lu and k\\ are the wave frequency and parallel wave number. Furthermore, 
the typical perpendicular velocity of the resonant electrons equals the electron thermal velocity, so that 
v± Vf <^ f l|. ThusEq. (26) is to be evaluated with U|| = ixi/k\\Vr and u±_ = 0. This gives 

#i = (27a) 

Pin U 

This efficiency is plotted in Fig. 7(a). Approximate fits for this function are given in Appendix B. 

On the other hand, for electron-cyclotron waves which interact through the Doppler- shifted cyclotron 
resonance oj — = Zf2, where Vl is the cyclotron frequency and / is the harmonic number, we have 
S = u_L. In this case we evaluate Eq. (26) at M|| — {lo ~ lVL)/k\\ Vr and m_l = to give 

P^ _ dWJdu-{llu^)dWJdn 

Pin" U ' ^^'"^ 

which is plotted in Fig. 7(b). 

Using Eqs. (27), it is possible to identify regions of high conversion efficiency of wave energy to 
electric energy, given the restriction of R small. Additionally, if the ohmic losses, y^/Psp^ small, 
then by Eq. (25), we see that the conversion of wave energy to poloidal field energy can be of high 
efficiency. This, in fact, is what has been achieved on the plt experiment, where conversion efficiencies 
of over 25% have been reported."^ 

An important practical consequence of the circuit equations derived here is that fast ramp-up rates, 
i.e., large 7, are possible at high density. In fact, these fast ramp-up rates are necessary for high energy 
conversion efficiencies at high density. This can be seen as follows: The efficiency, Poi/Pin> is a function 
of the dimensionless parameter wy, depending, in addition, only weakly on Z. For a given machine and 
a given wave phase velocity, the parameter u\\ depends only on the ratio, E/n; and the ramp-up rate, 
I, depends only on the dc electric field, E. Thus, the efficiency depends only on the ratio of I to n. It 
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has been observed experimentally on the PLT experiment that high efficiency of converting rf energy 
to magnetic field energy is possible at a low plasma density. Thus, we can predict that a similar high 
efficiency is possible in the event that the density and the ramp-up rate are scaled up together. In fact, for 
large ramp-up rates, high density can actually be desirable in that it impedes the production of runaways. 
Note that this window of desired density for a given ramp-up rate is counter to our intuition derived from 
steady-state considerations, where the larger the density the less the current-drive efficiency. 

There are several optimizations that one might wish to achieve in the ramp-up problem. One is 
to maximize the energy conversion efficiency, Pe\/P\n. A second is to minimize the ramp-up time, 
Tramp = I /I- The minimization of capital costs for the rf system, however, may demand that we 
minimize Prf , the rf power required to ramp-up a given current. 

We can express this more precisely with some convenient formulas. In the absence of the external 
source 14xt. the ramp-up rate may be written using Faraday's law as 

. 5£ MA 

(28) 



InRo/a s 

where a tokamak inductance L w i^qRo lni?o/a was assumed, and where E is the dc electric field in 
units of V/m. Note that the ramp-up rate depends hnearly on E and is almost independent of geometry 
(In i?o /aw 1). The amount of dissipated rf power required can than be written as 

(29) 

'/ ramp / m 

if we neglect both ohmic losses (F^/iisp) and the external source Pext- Thus, in extrapolating results 
to larger tokamaks (higher ^LP), we can maintain linear (in the required stored energy) power require- 
ments, with the same ramp-up time and the same efficiencies, if the density scales linearly with / and 
hence with /. Here, the wave phase velocities also remain the same, and in the event of the same temper- 
atures, the physics of the damping may be expected to be very similar, so that the percentage of incident 
rf power that is absorbed, 77, remains constant too. 

For example, using Fig. 7(a), we see that in PLT with n « 2 x 10^^ cm~'^, T w 1 keV, w/fcy « 
Qvt ~ \c, we find reported ramp-up rates of / = 120kA/s, 01 E 24mV/m. Here u\\ k, 1.4 and 
ideal efficiencies of about 33% at Z = 1 may be expected with little runaway production, consistent 
with experimental data. Also, consistent with the data would be somewhat higher Z, but then only if the 
confinement of runaways were not perfect. 

For reactor-grade tokamaks, say LP w 400 M J and 7 10 MA, a ramp-up time which is longer 
than that in PLT is desirable in order to minimize P^i and the capital cost of the rf system. For a 30 s 
ramp-up time, a density of 5 x 10^^ cm~^ renders the ratio I/n&s in the PLT experiment. Employing a 
similar spectrum of waves iuj/k\\ = jc) in a plasma of temperature also similar to the plt experiment 
(T « 1 keV) implies a similar rj (about 0.7). Thus, using Eq. (29), we see that Prf ~ 40 MW would be 
required. 

To summarize the tradeoffs here, we note that while Pel /Pin is minimized by considering only the 
ratio I/n, the minimization of Pif requires the Tramp be large. Thus, although very quick ramp-up 
rates are indeed achievable at high density, the capital costs for such a system are proportionately larger 
too. Balancing the desires for a quick ramp-up against those for low capital costs (low Prf) points to 
a parameter range of moderate density. Efficient ramp-up is only achieved when, in addition to the 
above restrictions, the temperature is moderate, since at high temperatures, V^/Rsp losses, neglected in 
Eq. (29), begin to dominate. The regime where these ohmic losses dominate may be identified by writing 

-^Sp ^ramp ^ramp 
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These losses represent only small corrections when y^/i?sp ^ -Prf, or, using Eq. (29), when 

For the reactor-grade example, Pei/Pin ~ 5, ^7 ~ 0.7, the above inequahty requires that the ramp-up 
time be longer than about j of the L/Rsp time. This restricts the temperature to somewhat less than 

2keV. 

Restricting the temperature during a period of intense rf injection (perhaps 40 MW) requires a small 
heat confinement time during the start-up operation. In the above example, this may be as small as 30 ms. 
Poor confinement during the start-up phase may be helpful from the standpoint of runaway buildup too. 
Even a small percentage (~ 1%) of reverse runaway s^^ can seriously impede ramp-up if the runaways are 
well confined. If the runaways are poorly confined, then higher percentages may be tolerated, allowing 
higher ramp-up rates and, consequently, higher energy conversion efficiencies for a given density. 

We are led thus to the following typical picture of rf ramp-up for pulsed tokamak operation. Start-up 
can proceed in a low density plasma^^ where the rf power is also used to initiate the plasma. Density 
and rf power, and the ramp-up rate, are increased concomitantly as the plasma is brought to interesting 
densities lO^'^-lO^'* cm^''. During this phase, the temperature is purposefully kept low, possibly through 
a deliberate degradation of the confinement of both runaway and thermal electrons. Hence the current is 
programmed to reach a large value prior to the density, and both reach large values prior to the tempera- 
ture. The final step, in which the reactor is brought to reactor-grade temperature, occurs after the current 
is ramped up and as a result of ceasing the deliberate degradation of confinement. 

VII. CONCLUSIONS 

In this paper we have written down a set of circuit equations that describe the dynamics of an rf- 
driven plasma. In arriving at these circuit equations, we systematically introduced approximations with 
a goal of characterizing the driven plasma by a small number of functions of few variables that retain 
the essential physics. Greater accuracy, possible at the price of more complex circuit equations, may 
be obtained as a natural extension of the development here. The identification and calculation here of a 
minimal set of transport functions, however, provided a suitable and manageable description for a large 
class of important problems. 

The calculations of the runaway function R and of the energy conversion function Ws together pin- 
point the preferred region for tokamak ramp-up operation. These functions depend only on the dimen- 
sionless parameter u. The separate contributions of runaway and stopped currents may be described 
using these functions of a single variable. The constitutive relations thus obtained are given by Eqs. (21). 
These equations are in a form both suitable for implementation in a transport code and amenable to 
obvious modification in the event that more complex runaway models are desired. 

There are several caveats to bear in mind in using these formulas. First, the time scale for variation of 
the dc electric fields has been assumed long compared to other scales of interest, such as the particle 
deceleration times. A violation of this scale separation would affect the normalizations through Vr. 
Second, knowledge of the rf spectrum is unlikely to be complete. This knowledge is necessary to give S, 
the rf-induced flux. Even if the incident rf energy is followed by ray-tracing codes, it remains possible 
that other waves may be present. These other waves might arise either due to asyrmnetries in the particle 
distribution functions or due to nonlinear effects associated with the incident spectrum. Third, particle 
transport across field lines was neglected in comparison to the effects along field lines. The neglect of 
these effects is possible for stopped electrons if they are stopped before they reach a flux surface with 
significantly different conditions {Vr different). For runaway electrons, these effects are always important 
in that they provide a model for the runaway loss. As discussed after Eqs. (21), such model may be 
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included through a natural modification of Eq. (21c). In the absence of one particularly compelling 
model, at present, for runaway loss we have left the modification of Eq. (21c) as an open issue. 

Finally, we should note that some of the most powerful conclusions of this paper occur in certain 
special cases. It is often the case that the rf spectrum is not only known, but also locahzed, which enables 
a particularly simple evaluation of the conversion efficiency, as in Eqs. (27). In the event of moderate 
electric fields, or spectra localized at moderate phase velocities, it may be that R = (no runaway 
production), and an accurate runaway loss model would not be needed. In the event that runaways are 
confined well, the spectrum must be chosen carefully to assure that R = 0. 
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Appendix A. LANGEVIN EQUATIONS 

Here we show that the conditional probability distribution g{u, r; u') for Eqs. (4) satisfies Eq. (3). 
The derivation follows those given in Refs. 18 and 25. Because the process described by Eqs. (4) is a 
Markoff process, g satisfies the Smolucowski equation^* 

g{u, T + At; u') = J d?v^' g{u, At; u")g{u", r; u') (Al) 

for all r > and At > 0. Let us define 

r(w, u, At) = ^(u + w, At; u). 
Subtracting ^'(u, r; u') from Eq. (Al) gives 

gin, T + At; u') - ff(u, r; u') = j d'u" [r(u - u", u", AT)g(u", r; u') 

-r(u-u",u,AT)ff(u,T;u')], (A2) 

where, because of the normalization condition for probabilities, the second term in the integral may be 
reduced to g{u, r; u'). If we change the variable of integration to w = u — u", the right-hand side of 
Eq. (A2) becomes 

J d'^w [r(w, u — w, AT)sf(u — w, r; u') — r(w, u, Ar)(?(u, r; u')] . 

For small At, the function r(w, u, At) is highly localized about w = 0. We may therefore expand the 
first term in the integral, assuming that w is much smaller than u, to give 

d 

r(w, u - w, AT)fif(u - w, t; u') w r(w, u, AT)fif(u, t; u') - w • -g^r{w, u, AT)fif(u, t; u') 

1 

2 ouou 

Using this approximation in Eq. (A2), integrating by parts, dividing by At and taking the hmit At 0, 
we find as the equation for ^(u, t; u') 

d d 
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where 

A(u)= lim 



At^o At 
B(u) = lim ^— 

^ ' Ar^O 2At 

and 

(Au) = y" wr(w,u, Ar)d^w, 



( AuAu) = y WW r(w, u, At) d^w. 



Thus (Au) is the average value of u(t + Ar) — u(r) given that u(r) = u (and similarly for (AuAu)). 
These quantities may be calculated directiy from Eqs. (4) assuming that At is sufficiently small that u 
does not change appreciably. We then obtain 

<.r+AT 

(A«) 



i+M)AT, 

r+Ar / 1 _ ,,f^/\2 ' 



l + Z 1 - , ^ 
^ At, 

r+Ar pT+Ar 



(A/iA/x)=y"" "(A(T')^(r")) rfr'dT" + 0(AT2) 

«^(1-A^^)At, 
(AuAu) = (AuA^) = (A/iAu) = ©(At^). 

Here we have made use of the properties of A given in Eqs. (5). Writing Eq. (A3) in spherical coordinates 
and substituting for the non-zero components of A and B, we obtain 

d Id 1 + Z d 2, d 

= d^^' + ^d^' + ^djl^' - ^ ^T^'- ^^^^ 

This is the same equation as Eq. (3). Furthermore, from the definition of g as a conditional probabihty, 
the initial condition for Eq. (A4) is also the same as for Eq. (3), namely, gin, t; u') = (5(u — u'). 



Appendix B. NUMERICAL FITS 



In this appendix, we give approximations for some of the important functions we have calculated. 
These are suitable for incorporating into modeUng codes. The approximations were found by choosing 
a suitable analytic form containing several undetermined coefficients and adjusting those coefficients in 
order to minimize the maximum relative error. The technique for carrying out this procedure is described 
in Hastings' classic work.^^ The fits were made to the numerical data presented in Sec. IV. This data 
itself contains errors due to the numerical methods used. The main source of error is due to the finite size 
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of the numerical mesh and it is estimated that this introduces errors on the order of a percent. However, 
near u = 0, the relative error in the numerical data for Wg and its derivative becomes large because 
Ws = 0{u'^). Thus for u < 0.5, the fits were made using the following analytical approximation instead 
of the numerical data: 

^_ (2 + Z + 3/x2)w6 ^ 2((24+19Z + 3Z2)/x + (9 + Z)/x3)u8 



Z + 5 3(3 + Z)(5 + Z) {S + Z){5 + Z){7 + 3Z){9 + Z) 

/ 1041 + 1864Z+11892^+316Z'''+30Z* \ lO 

_ V +10(417+497Z+181Z2+21Z3)^2_,_5(g+2)(13^32)^4 

~ 5(2 + Z){3 + Z){5 + Z){7 + 3Z){9 + Z){13 + 3Z) ' 

This result was obtained by solving Eq. (12) for small u using MACSYMA.^'' (The first two terms in this 
expansion are those derived by Fisch.^^) 

For each of the functions approximated, we give the analytic form of the approximation, the range 
in which it is valid, a table of coefficients, and the maximum relative error. The approximations should 
not be used outside the range given. Also, note that the relative error quoted is the error in fitting the 
approximation to the numerical data which itself is in error by about a percent. 

For /U = 1 and 1.4 < u < 8, the runaway probability R is approximated by 

i?(u,M=l) = expf ^-=°"'("-^)^ ^ 



where &i = 1 and the other coefficients o, and bi are given in Table 1. The maximum relative error is 1%. 
For /u = 1 and 1 < m < 1.4, the same approximation may be used with small absolute error but large 
relative error. For u < 1 and all /x we have R = identically. 

For IJ. = 1 and < u < 5, the energy imparted to the electric field Wg by stopped electrons is 
approximated by 



Wg{u,^i = l) 



where bo = I and the other coefficients a, and bi are given in Table II. The maximium relative error is 
2%. For n = —1 and < m < 1, Wg is approximated by 



Ws{u,n= -1) = ^tti 



,2i 



where the coefficients are given in Table III and the maximium relative error is 1.5%. 
For = 1 and < u < 5, the function {dWs / du) /u is is approximated by 

where bo = 1 and the other coefficients and bi are given in Table IV. The maximium relative error is 
5%. For = — 1 and < u < 1, {dWs/du)/u is approximated by 

1 d 

--W^,(t.,M=-l) = E«'^"'* 



u du 

2—1 



where the coefficients aj are given in Table V and the maximium relative error is 3%. 
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Tables 

TABLE I. Coefficients for approximation to R{u, /i = 1). 



z 


ao 


ai 


02 


03 


b2 


bs 


1 


-3.68063 


4.23913 


-4.55894 


-0.39755 


-1.22774 


1.41450 


2 


-4.97636 


-16.09015 


0.83188 


0.21737 


6.84615 


-0.98649 


5 


-4.27687 


-4.33629 


0.30338 


0.05697 


3.21315 


-0.47749 


10 


-4.94597 


-1.53482 


0.10112 


0.03087 


2.45288 


-0.36896 



TABLE II. Coefficients for approximation to Ws{u, fi — 1). 



Z_ 02 03 £4 bi 62 63 

1 0.16612 -0.01495 0.00775 0.37136 0.02240 0.01645 
2 0.14200 -0.04048 0.01145 0.12253 0.00384 0.02440 
5 0.09880 -0.05152 0.01113 -0.19484 0.00559 0.02362 
10 0.06537 -0.03895 0.00738 -0.32456 0.02797 0.01526 



TABLE III. Coefficients for approximation to Ws{u, = —1). 



Z 02 03 04 05 

1 -0.16483 -0.13420 0.15346 -0.24314 

2 -0.14186 -0.09297 0.06661 -0.12870 

5 -0.09975 -0.04781 0.00606 -0.03545 

10 -0.06651 -0.02797 -0.00247 -0.00934 



TABLE IV. Coefficients for approximation to {dWs/ du) / u fox ijl = 1. 



z 


Ol 


02 


03 


bi 


b2 


bi 


1 


0.66445 


-0.36032 


0.07328 


0.17769 


-0.25452 


0.07278 


2 


0.56760 


-0.38984 


0.08634 


-0.04019 


-0.24673 


0.08508 


5 


0.39906 


-0.32879 


0.07670 


-0.28281 


-0.16275 


0.07436 


10 


0.27028 


-0.23261 


0.05272 


-0.39140 


-0.07526 


0.04981 



TABLE V. Coefficients for approximation io {dWs/du)/ u for ji = — 1 . 



Z ai a2 03 04 

1 -0.63673 -1.39960 3.37662 -4.23684 

2 -0.55777 -0.80763 1.43144 -2.03866 

5 -0.39704 -0.33811 0.23607 -0.51011 

10 -0.2()G00 -0.17342 0.0189(i -0.13349 
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FIG. 1. The runaway probability -R(u) foxZ= 1. Parts (a) and (b) show R on two different scales. In (a) 
the contours are equally spaced at intervals of 0.05. In (b) the lowest 7 contours are geometrically spaced 
at intervals of 10^/^ between 10~^ and 10~^; the remaining contours are equally spaced at intervals of 
0.05 as in (a). 
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FIG. 2. i?(u|| , u_L = 0) for Z=l,2, 5, and 10. 
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FIG. 3. The current j(u, r) for u = 5u|| and Z = 1; (a) the total current j; (b) the stopped current jg; (c) 
the runaway current jV- For u = 5u||, approximately 32% of the electrons run away. 
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FIG. 4. The energy imparted to the electric field by the stopped particles Ws (u) for Z = 1. The innermost 
contours are equally spaced at intervals of 0.005 between —0.05 and 0.05. The remaining contours are 
equally spaced at intervals of 0.05. Part (a) shows the results of numerically solving Eq. (12); part (b) 
shows Ws from the hot-conductivity theory Eq. (13). 
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FIG. 6. The perpendicular energy of the runaways f_Lr(u) for Z = 1. The contours are equally spaced 
at intervals of 0.5. 
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